#!/bin/tcsh -ef
####
#
#############################################################################
#                                                                           #
# Title: Xstal Filter (Experimental)                                        #
#                                                                           #
# (C) 2dx.org, GNU Plublic License.                                         #
#                                                                           #
# Created..........: 04/24/2019                                             #
# Last Modification: 04/24/2019                                             #
# Author...........: focus-em.org                                           #
#                                                                           #
#############################################################################
#
# SORTORDER: 37
#
# MANUAL: This script determines the <B>reciprocal crystal lattice</B> in a Fourier transformation of the image. Two different algorithms are offered. The first algorithm <I>FindLattice</I> makes use of a-prior knowledge of the real-space dimensions of the crystal lattice, and of the tilt geometry. The second algorithm <I>GetLattice</I> performs a peak list difference calculation and determines the most likely lattice from the power-spectrum.
# 
# MANUAL: <B><U>FindLattice</U></B>
#
# MANUAL: This algorithm determines the lattice, using the a-priori knowledge of the expected lattice and the approximate tilt geometry (from the defocus determination). Once the lattice is determined, its distortion is used to determine the tilt geometry again, using the MRC program EMTILT (D. Agard). With tilt angles higher than 25 degrees, the lattice distortion will usually allow a better determination of the tilt geometry than from a defocus gradient. 
#
# MANUAL: <B><U>GetLattice</U></B>
#
# MANUAL: This algorithm <I>Difference Lattice</I> does not require a priori knowledge of the real-space lattice. 
#
# MANUAL: <i>2dx_diffLattice</i> takes advantage of the shifted origin peak list. The results of this shifting are identical to those produced by weighted, difference-vector methods --with the added benefit of using an adaptive criteria for peak-merging as opposed to standard Gaussian determination.
#
# MANUAL: <i>2dx_diffLattice</i> begins by finding all peaks with radii less than the maximum radius of the 8 strongest peaks from the list. It follows from the basis for the shifted origin that the strongest peak will necessarily fall on the stronger of the two lattices. Every possible vector pair is then used to generate a lattice which is then compared against the available peak list. Lattice error is computed by transforming peaks into Miller indices, then measuring their divergence from integer values, producing a result which is then multiplied by the strength of the point in question. For lattice pairs with an angle of 10 degrees or less between them, a further weighting factor of Exp(10-Theta) is multiplied, with theta the measured angle between the vectors in degrees. (This somewhat arbitrary weighting is used to prevent extremely high density lattices from overpowering sparse, yet more accurate, ones.)
#
# MANUAL: The vector pair which generates the lowest error lattice is then used as a Miller basis. Using this basis, a Miller index is then assigned to any peak which lies within epsilon=0.15*sqrt(2) of a point on the lattice. A peak assigned to an index is replaced if another peak is found to lie closer to the lattice point. Peaks which are not assigned lattice points are discarded from further consideration.
#
# MANUAL: Using the generated Miller index - peak position pairs, a least squares fit is then performed to determine a candidate lattice which is then used as the basis for iterative refinement, according to the above procedure, to produce the <B>final lattice</B>.
#
# MANUAL: <HR>
#
#
# DISPLAY: movie_stackname
# DISPLAY: lattice
# DISPLAY: backuplattice
# DISPLAY: defocus
# DISPLAY: defocusfield
# DISPLAY: RESMAX
# DISPLAY: realang
# DISPLAY: realcell
# DISPLAY: SYM
# DISPLAY: tempkeep
# DISPLAY: delta
# DISPLAY: streakfactor
# DISPLAY: FirstPeakNum
# DISPLAY: peakNum
# DISPLAY: testHand
# DISPLAY: InnerExclusionRadius
# DISPLAY: comment
# DISPLAY: det_lattice
# DISPLAY: det_tilt_lat
# DISPLAY: do_lattice_algorithm
# DISPLAY: laterror_abserror
# DISPLAY: laterror_rmsderror
# DISPLAY: laterror_peaksused
# DISPLAY: laterror_nodepeakdensity
# DISPLAY: invmask_doit
# DISPLAY: invmask_radius
# DISPLAY: invmask_type
# DISPLAY: invmask_spots
# DISPLAY: stepdigitizer
# DISPLAY: IQMAX
# DISPLAY: TLTAXIS
# DISPLAY: TLTANG
# DISPLAY: spotlist_RESMAX
# DISPLAY: ALAT
#
#$end_local_vars
#
#
set bin_2dx = ""
set proc_2dx = ""
#
set SYN_Unbending = ""
set PHASEORI_done = ""
#
set movie_stackname = ""
set imagenumber = ""
set realcell = ""
set lattice = ""
set TLTAXIS = ""
set TLTANG = ""
set TLTAXA = ""
set TAXA = ""
set TANGL = ""
set imagesidelength = ""
set magnification = ""
set stepdigitizer = ""
set tempkeep = ""
set revhk = ""
set realang = ""
set RESMIN = ""
set RESMAX = ""
set ALAT = "" 
set det_tilt_lat = ""
set delta = ""
set streakfactor = ""
set FirstPeakNum = ""
set peakNum = ""
set testHand = ""
set regenPL = ""
set InnerExclusionRadius = ""
set det_lattice = ""
set do_lattice_algorithm = ""
set laterror_abserror = ""
set laterror_rmsderror = ""
set laterror_peaksused = ""
set laterror_nodepeakdensity = ""
set MASKING_done = ""
set LATTICE_done = ""
set defocus = ""
set invmask_doit = ""
set invmask_radius = ""
set invmask_type = ""
set invmask_spots = ""
set stepdigitizer = ""
set IQMAX = ""
set spotlist_RESMAX = ""
#
#$end_vars
#
set scriptname = process_XFilter
#
\rm -f LOGS/${scriptname}.results
#
if ( ${invmask_doit} == "n" ) then
  protest "Not running."
endif
#  
set basel = `uname -a | grep ethz.ch | wc -l`
if ( ${basel} == '0' ) then
#  protest "ERROR: In preparation, do not use."
endif
#
source ${proc_2dx}/initialize
#
echo "<<@evaluate>>"
#
set docfile = 'peaks_xy_final.dat'
set docfileraw = 'peaks_xy.dat'
set iradius = '8'
set imaxrad = '20'
#
set date = `date`
echo date = ${date}
#
if ( ${det_lattice}x == 'x' ) then
  set det_lattice = 'y'
  echo "set det_lattice = ${det_lattice}" >> LOGS/${scriptname}.results
  echo ":: det_lattice corrected to ${det_lattice}"
endif
#
if ( x${TLTAXIS} == "x-" || x${TLTAXIS} == "x--" ) then
  set TLTAXIS = "0.0"
endif
if ( x${TLTANG} == "x-" || x${TLTANG} == "x--" ) then
  set TLTANG = "0.0"
endif
#
set input_image = ${movie_stackname}
#
if ( ! -e ${input_image}.mrc ) then
  ${proc_2dx}/protest "ERROR: ${input_image}.mrc does not exist."
endif
#
if ( -e ${movie_stackname}_Sum.mrc ) then
  set input_image = ${movie_stackname}_Sum
  echo "# IMAGE-IMPORTANT: ${movie_stackname}_Sum.mrc <DriftCor image (2D, no DW)>" >> LOGS/${scriptname}.results
  echo "# IMAGE-IMPORTANT: ${movie_stackname}_Sum_fft.mrc <DriftCor image FFT (2D, no DW)>" >> LOGS/${scriptname}.results
  echo "# IMAGE: ${movie_stackname}.mrc <DriftCor image (2D, with DW)>" >> LOGS/${scriptname}.results
  echo "# IMAGE: ${movie_stackname}_fft.mrc <DriftCor image FFT (2D, with DW)>" >> LOGS/${scriptname}.results
else
  echo "# IMAGE-IMPORTANT: ${movie_stackname}.mrc <DriftCor image (2D, with DW)>" >> LOGS/${scriptname}.results
  echo "# IMAGE-IMPORTANT: ${movie_stackname}_fft.mrc <DriftCor image FFT (2D, with DW)>" >> LOGS/${scriptname}.results
endif
#
setenv IN ${movie_stackname}_Sum.mrc
${bin_2dx}/2dx_header.exe
set imagesidelength = `cat 2dx_header.out | head -n 4 | tail -n 1 | cut -c43-52`
set imagesidelengthx = `cat 2dx_header.out | head -n 4 | tail -n 1 | cut -c43-52`
set imagesidelengthy = `cat 2dx_header.out | head -n 4 | tail -n 1 | cut -c53-62`
#
echo "::Defocus = "${defocus}
#
if ( ${det_lattice} == 'n' && ${det_tilt_lat} == 'n' ) then
  ${proc_2dx}/linblock "Skipping. Not determining lattice."
  exit
endif
#
\rm -f TMP456789.tmp
#
echo "<<@progress: 1>>"
#
if ( ${det_lattice} == 'y' ) then
  #
  set ialg = `echo ${do_lattice_algorithm} | cut -c1-1` 
  #
  if ( ${ialg} == "0" ) then
    #
    #############################################################################
      ${proc_2dx}/linblock "Running FindLattice"
    #############################################################################
    #
    if ( "${realcell}" == "100.0,100.0" || "${realcell}" == "100,100" ) then
      echo ":: "
      echo ":: "
      echo "::                        WARNING: You are using FindLattice,"
      echo "::                     but did you define a real lattice already?" 
      echo ":: "
      echo "::        (The current value of ${realcell} looks like you didn't do that yet.)"
      echo "::       You might need to first determine the real unit cell length and angle,"
      echo "::                            or use GetLattice instead."
      echo ":: "
      echo ":: "
      echo "#WARNING: Warning: For FindLat you need a valid real lattice first." >> LOGS/${scriptname}.results
    endif
    #
    if ( "${realcell}" == "0,0" || "${realcell}" == "0.0,0.0" ) then
      echo ":: "
      echo ":: "
      echo "::                        ERROR: You are using FindLattice,"
      echo "::                     but there is no real lattice defined." 
      echo ":: "
      echo "::            (The current value for the real lattice is ${realcell}.)"
      echo "::        You need to first determine the real unit cell length and angle."
      echo "::      You can do that manually, or use GetLattice instead. After that, run"
      echo "::     the special script Evaluate Lattice to determine the real cell values."
      echo ":: "
      echo ":: "
      echo "#WARNING: ERROR: For FindLat you need a valid real lattice first." >> LOGS/${scriptname}.results
      ${proc_2dx}/protest "Aborting."
    endif
    #
  else
    #
    #############################################################################
    ${proc_2dx}/linblock "Running GetLattice"
    #############################################################################
    #
  endif
  echo "<<@progress: 5>>"
  #
  if ( ${regenPL} == "y" || ( ! -e ${docfile} ) || ( ! -e ${docfileraw} ) ) then
    #
    #############################################################################
    ${proc_2dx}/linblock "2dx_taperedge: To taper the edges to prevent stripes in the FFT"
    #############################################################################
    #
    set largeimage = `echo ${imagesidelength} | awk '{ if ( $1 > 6000 ) { s = 1 } else { s = 0 } } END { print s }'`
    #
    setenv IN ${movie_stackname}_Sum.mrc
    \rm -f     SCRATCH/${movie_stackname}_Sum_taper.mrc
    setenv OUT SCRATCH/${movie_stackname}_Sum_taper.mrc
    ${bin_2dx}/2dx_taperedgek.exe << eot
10,10,100,10       ! IAVER,ISMOOTH,ITAPER 
eot
    #
    echo "<<@progress: 10>>"
    #
    #############################################################################
    ${proc_2dx}/linblock "2dx_peaksearch: To generate a completed peak list ..."
    #############################################################################
    #
    echo "<<@progress: 20>>"
    #
    \rm -f ${docfile}
    \rm -f ${docfileraw}
    \rm -f 2dx_peaksearch-amp_before_LHpass.mrc
    \rm -f 2dx_peaksearch-amp_LHpass.mrc
    #
    echo "FirstPeakNum = ${FirstPeakNum}"
    echo "peakNum = ${peakNum}"
    echo "InnerExclusionRadius = ${InnerExclusionRadius}"
    echo "imagesidelength = ${imagesidelength}"
    echo "streakfactor = ${streakfactor}"
    ${bin_2dx}/2dx_peaksearch.exe SCRATCH/${movie_stackname}_Sum_taper.mrc ${FirstPeakNum} ${peakNum} ${InnerExclusionRadius} ${imagesidelength} ${streakfactor}
    #
    if ( ${tempkeep} == 'y' ) then
      #
      \mv -f 2dx_peaksearch-non-masked_image.mrc SCRATCH/pks-non-masked_image.mrc
      if ( -e 2dx_peaksearch-amp_before_LHpass.mrc ) then
        \mv -f 2dx_peaksearch-amp_before_LHpass.mrc SCRATCH/pks-amp_before_LHpass.mrc
      endif
      if ( -e 2dx_peaksearch-amp_LHpass.mrc ) then
        \mv -f 2dx_peaksearch-amp_LHpass.mrc SCRATCH/pks-amp_LHpass.mrc
      endif
      \mv -f 2dx_peaksearch-masked_image.mrc SCRATCH/pks-masked_image.mrc
      \rm -f 2dx_peaksearch-phase.mrc
      echo "# IMAGE: SCRATCH/${movie_stackname}_Sum_taper.mrc <Edge-Tapered Image>" >> LOGS/${scriptname}.results
      echo "# IMAGE: SCRATCH/pks-non-masked_image.mrc <Powerspectrum>" >> LOGS/${scriptname}.results
      echo "# IMAGE: SCRATCH/pks-amp_before_LHpass.mrc <Pre Band-Pass Powerspectrum>" >> LOGS/${scriptname}.results
      echo "# IMAGE: SCRATCH/pks-amp_LHpass.mrc <Post Band-Pass Powerspectrum>" >> LOGS/${scriptname}.results
      echo "# IMAGE: SCRATCH/pks-masked_image.mrc <Masked and Band-Passed Powerspectrum>" >> LOGS/${scriptname}.results
      #
    else
      \rm -f SCRATCH/${movie_stackname}_Sum_taper.mrc
      \rm -f 2dx_peaksearch-non-masked_image.mrc
      \rm -f 2dx_peaksearch-amp_LHpass.mrc
      \rm -f 2dx_peaksearch-masked_image.mrc
      \rm -f 2dx_peaksearch-phase.mrc
    endif
    #
    \mv -f 2dx_peaksearch-average.mrc SCRATCH/pks-average.mrc
    echo "# IMAGE-IMPORTANT: SCRATCH/pks-average.mrc <Origin-Shifted Average Powerspectrum>"  >> LOGS/${scriptname}.results
    #
  else
    #############################################################################
    ${proc_2dx}/linblock "2dx_peaksearch: Skipping"
    #############################################################################
    echo "# IMAGE: SCRATCH/${movie_stackname}_Sum_taper.mrc <Edge-Tapered Image>" >> LOGS/${scriptname}.results
    echo "# IMAGE: SCRATCH/pks-non-masked_image.mrc <Powerspectrum>" >> LOGS/${scriptname}.results
    echo "# IMAGE: SCRATCH/pks-amp_LHpass.mrc <Band-Passed Powerspectrum>" >> LOGS/${scriptname}.results
    echo "# IMAGE: SCRATCH/pks-masked_image.mrc <Masked and Band-Passed Powerspectrum>" >> LOGS/${scriptname}.results
    echo "# IMAGE-IMPORTANT: SCRATCH/pks-average.mrc <Origin-Shifted Average Powerspectrum>" >> LOGS/${scriptname}.results
    #
  endif
  #
  echo "<<@evaluate>>"
  #
  if ( ${ialg} == "0" ) then
    #
    echo "<<@progress: 35>>"
    #
    #############################################################################
    ${proc_2dx}/linblock "2dx_findlat: To determine a lattice ..."
    #############################################################################
    #
    \rm -f 2dx_findlat.out
    # 
    set tltAx = ${TLTAXIS}
    set tltAn = ${TLTANG}
    #
    if ( ${testHand} == 'n' ) then
      set itesthand = 0
    else
      set itesthand = 1
    endif
    #  
    echo "starting findlat with:"
    echo "${realcell}"
    echo "${realang}"
    echo "${tltAx}"
    echo "${tltAn}"
    echo "${imagesidelength}"
    echo "${magnification}"
    echo "${stepdigitizer}"
    echo "${docfile}"
    echo "0.1,0.008,8,${delta},16,300,${itesthand}"
    echo " "
    #
  ${bin_2dx}/2dx_findlat.exe << eot
${realcell}
${realang}
${tltAx}
${tltAn}
${imagesidelength}
${magnification}
${stepdigitizer}
${docfile}
0.1,0.008,8,${delta},16,300,${itesthand} !AngleStepSize,MagVariation,Nmin,Delta,MaxHK,MaxSpot
eot
    #
    echo "findlat finished."
    #
    if ( ${tempkeep} == 'n' ) then
      #\rm -f peaks_xy_final.dat
    endif
    #
    if ( ! -e 2dx_findlat.out ) then
      ${proc_2dx}/protest "2dx_findlat: ERROR occured."
    endif
    #
    set tmplat = `cat 2dx_findlat.out | head -n 1`
    set firstscore  = `cat 2dx_findlat.out | head -n 2 | tail -n 1`
    set newlat = `echo ${tmplat} | sed 's/ /,/g'`
    set oldlat = `echo ${lattice}`
    #
    echo "set lattice = "\"${newlat}\" >> LOGS/${scriptname}.results
    # "
    set lattice = ${newlat}
    echo "::  Old lattice is ${oldlat}"
    echo "::  New lattice is ${lattice} (Score = ${firstscore})"
    #
    echo " " >> History.dat
    echo ":Date: ${date}" >> History.dat
    echo "::Lattice is ${lattice}" >> History.dat
    #
  else
    #
    #############################################################################
    ${proc_2dx}/linblock "2dx_getlat: To determine a lattice ..."
    #############################################################################
    #
    set tltAx = ${TLTAXIS}
    set tltAn = ${TLTANG}
    #
    set lattice=`${bin_2dx}/2dx_getlat.exe ${docfile} | tail -n 1`
    #
    echo "set lattice = "\"${lattice}\" >> LOGS/${scriptname}.results
    # "
    echo " " >> History.dat
    echo ":Date: ${date}" >> History.dat
    echo "::Lattice is ${lattice}" >> History.dat
    #
    if ( ${tempkeep} == 'n' ) then
      #\rm -f peaks_xy_final.dat
    endif
    #
    echo "::  New lattice is ${lattice}"
  endif
  #
  echo "<<@progress: 50>>"
  echo "<<@evaluate>>"
  #
  #############################################################################
  ${proc_2dx}/linblock "2dx_mmboxa: To get the spot list from the FFT of the image."
  #############################################################################
  #
  set IQMAX_local = ${IQMAX}
  set imagecenterx = `echo ${imagesidelengthx} | awk '{ s = int( $1 / 2 ) } END { print s }'`
  set imagecentery = `echo ${imagesidelengthy} | awk '{ s = int( $1 / 2 ) } END { print s }'`
  #
  # Cutoff radius is:  
  #
  #                       Size * stepsize * 10000
  # Cutoff radius = ---------------------------------
  #                 Resolution[A] * Magnification
  #
  #
  # Example:
  #
  #    1024 * 17.0 * 10000
  #  --------------------- = 241
  #   16 * 45000
  #
  #
  #
  set RINNER = 1
  set ROUTER = `echo ${imagesidelength} ${stepdigitizer} ${magnification} ${spotlist_RESMAX} | awk '{ s = ( $1 * $2 *10000 ) / ( $3 * $4 ) } END { print s }'`
  ${proc_2dx}/linblock "Outer radius = ${ROUTER}"
  #
  echo RINNER ROUTER 0 realcell ALAT realang ! RINNER,ROUTER,IRAD,A,B,W,ABANG
  echo ${RINNER} ${ROUTER} 0 ${realcell} ${ALAT} ${realang} ! RINNER,ROUTER,IRAD,A,B,W,ABANG
  #
  #
  \rm -f ${movie_stackname}.spt
  ${bin_2dx}/2dx_mmboxa.exe << eot
${movie_stackname}_fft.mrc
${imagenumber}
Y                               ! Use grid units?
Y                               ! Generate grid from lattice?
N                               ! Generate points from lattice?
2 3 0 50 50 12 12               ! IPIXEL,IOUT,NUMSPOT,NOH,NOK,NHOR,NVERT
${movie_stackname}.spt
${IQMAX_local}                              ! IQMAX
${imagecenterx} ${imagecentery}         ! XORIG,YORIG
${RINNER} ${ROUTER} 0 ${realcell} ${ALAT} ${realang} ! RINNER,ROUTER,IRAD,A,B,W,ABANG
${lattice}                         ! Lattice vectors
eot
  #
  #
    #############################################################################
    #############################################################################
    #############################################################################
    #############################################################################
else
  if ( "${lattice}" == "0.0,100.0,100.0,0.0" || "${lattice}" == "0.0,0.0,0.0,0.0" ) then
    ${proc_2dx}/protest "ERROR: you need to determine the reciprocal lattice first."
  endif
  if ( "${realcell}" == "100.0,100.0" || "${realcell}" == "100,100" ) then
    echo "::        (The current value of ${realcell} looks like you did not do that yet.)"
    echo "::       You might need to first determine the real unit cell length and angle,"
    ${proc_2dx}/protest "ERROR: Aborting."
  endif
  if ( "${realcell}" == "0,0" || "${realcell}" == "0.0,0.0" ) then
    echo "::            (The current value for the real lattice is ${realcell}.)"
    echo "::        You need to first determine the real unit cell length and angle."
    echo "::      You can do that manually, or use GetLattice instead. After that, run"
    echo "::     the special script Evaluate Lattice to determine the real cell values."
    ${proc_2dx}/protest "Aborting."
  endif
endif
#############################################################################
#############################################################################
#############################################################################
#############################################################################
#
#
#############################################################################
${proc_2dx}/linblock "2dx_laterror to find fitness of current lattice."
#############################################################################
#
${bin_2dx}/2dx_laterror.exe "${lattice}" ${docfile} > SCRATCH/2dx_laterror.tmp
cat SCRATCH/2dx_laterror.tmp
set laterror_abserror=`cat SCRATCH/2dx_laterror.tmp | tail -n 7 | head -n 1`
set laterror_rmsderror=`cat SCRATCH/2dx_laterror.tmp | tail -n 5 | head -n 1`
set laterror_peaksused=`cat SCRATCH/2dx_laterror.tmp | tail -n 3 | head -n 1 | cut -f1 -d " "`
set laterror_nodepeakdensity=`cat SCRATCH/2dx_laterror.tmp | tail -n 1`
#
echo ": "
echo "::Absolute Error: ${laterror_abserror}"
echo "::RMSDn ERROR: ${laterror_rmsderror}"
echo "::Peaks Used: ${laterror_peaksused}"
echo "::Nodes/Peak Density: ${laterror_nodepeakdensity}"
#
echo "set laterror_abserror = ${laterror_abserror}" >> LOGS/${scriptname}.results
echo "set laterror_rmsderror = ${laterror_rmsderror}" >> LOGS/${scriptname}.results
echo "set laterror_peaksused = ${laterror_peaksused}" >> LOGS/${scriptname}.results
echo "set laterror_nodepeakdensity = ${laterror_nodepeakdensity}" >> LOGS/${scriptname}.results
#
#############################################################################
${proc_2dx}/linblock "LABEL - to create flat value=1 FFT file"
#############################################################################
#
\rm -f SCRATCH/TMPvalue1_fft.mrc
#
${bin_2dx}/labelh.exe << eot
${movie_stackname}_fft.mrc
9               ! Fill FFT with Amplitude
SCRATCH/TMPvalue1_fft.mrc
1.0
eot
  #
  echo "# IMAGE: SCRATCH/TMPvalue1_fft.mrc <Value 1 FFT>" >> LOGS/${scriptname}.results
  echo "<<@progress: +5>>"
  #
  #
  #############################################################################
  ${proc_2dx}/linblock "MASKTRAN - to delete spots of lattice"
  #############################################################################
  #
  setenv IN  SCRATCH/TMPvalue1_fft.mrc
  setenv OUT SCRATCH/TMPvalue1_msk_fft.mrc
  setenv SPOTS ${movie_stackname}.spt
  #
  \rm -f SCRATCH/TMPvalue1_msk_fft.mrc
  set ishape = 1
  set radius = 5
  set rmax = 11000
  echo rmax = ${rmax}
  if ( ${invmask_type} == "0" ) then
    # mask with circular mask
    set ishape = 1
    set radius = ${invmask_radius}
  endif
  if ( ${invmask_type} == "1" ) then
    # mask with Gaussian mask
    set ishape = 2
    set radius = ${invmask_radius}
  endif
  if ( ${invmask_type} == "2" ) then
    # mask with rectangular mask
    set ishape = 3
    set radius = "${invmask_radius},${invmask_radius}"
  endif
  #
  set itype = 1
  #  0=From Spotlist;1=All Spots
  if ( ${invmask_spots} == "0" ) then
    set itype = 1
  endif
  if ( ${invmask_spots} == "1" ) then
    set itype = 0
  endif
  #
  echo ":: Masking with first lattice: Now launching masktrana with ishape=${ishape}:" 
  cat << eof
::  ${bin_2dx}/2dx_masktrana.exe << eot
:: ${ishape} F T F ! ISHAPE=1(CIRC),2(GAUSCIR),3(RECT)HOLE,IAMPLIMIT(T or F),ISPOT,IFIL
:: ${radius}       ! RADIUS OF HOLE IF CIRCULAR, X,Y HALF-EDGE-LENGTHS IF RECT.
:: ${lattice} -50 50 -50 50 ${rmax} ${itype} !A/B X/Y,IH/IK MN/MX,RMAX,ITYPE
:: eot
eof
  #
  ${bin_2dx}/2dx_masktrana.exe << eot
${ishape} F T F ! ISHAPE=1(CIRC),2(GAUSCIR),3(RECT)HOLE,IAMPLIMIT(T or F),ISPOT,IFIL
${radius}       ! RADIUS OF HOLE IF CIRCULAR, X,Y HALF-EDGE-LENGTHS IF RECT.
${lattice} -50 50 -50 50 ${rmax} ${itype} !A/B X/Y,IH/IK MN/MX,RMAX,ITYPE
eot
  #
  echo "# IMAGE: SCRATCH/TMPvalue1_msk_fft.mrc <Masked Value 1 FFT>" >> LOGS/${scriptname}.results
  #
  echo "<<@progress: +5>>"
  #
  #############################################################################
  ${proc_2dx}/linblock "LABEL - to invert masked FFT pattern"
  #############################################################################
  #
  \rm -f SCRATCH/TMPvalue2_fft.mrc
  #
  ${bin_2dx}/labelh.exe << eot
SCRATCH/TMPvalue1_msk_fft.mrc
2               ! Linear OD stretch  ( y = mx + b )
SCRATCH/TMPvalue2_fft.mrc
-1.0 1.0
0
eot
  #
  echo "# IMAGE: SCRATCH/TMPvalue2_fft.mrc <Inverted masked Value 1 FFT>" >> LOGS/${scriptname}.results
  echo "<<@progress: +5>>"
  #
  #############################################################################
  ${proc_2dx}/linblock "TWOFILE - to multiply masked pattern with FFT of image"
  #############################################################################
  #
  setenv IN1 ${movie_stackname}_fft.mrc
  setenv IN2 SCRATCH/TMPvalue2_fft.mrc
  \rm -f     SCRATCH/${movie_stackname}_invmask_fft.mrc
  setenv OUT SCRATCH/${movie_stackname}_invmask_fft.mrc
  ${bin_2dx}/2dx_twofile.exe << eot
1               ! ICOMB (multiply together)
2 0 0 ${imagecenterx} ${imagecentery} ! IORIGIN,OXA,OYA,OXB,OYB  Origin shifts to FFTs
eot
  #
  echo "<<@progress: +5>>"
  #
  echo "# IMAGE-IMPORTANT: SCRATCH/${movie_stackname}_invmask_fft.mrc <Inverse Fourier-filtered Image FFT (2D, with DW)>" >> LOGS/${scriptname}.results
  #
  #############################################################################
  #                                                                           #
  ${proc_2dx}/linblock "clip fft - to calculate Fourier filtered image"
  #                                                                           #
  #############################################################################
  #
  \rm -f                                                              ${movie_stackname}_invflt.mrc
  ${dir_imod}/bin/clip fft SCRATCH/${movie_stackname}_invmask_fft.mrc ${movie_stackname}_invflt.mrc
  #
  # ${dir_imod}/bin/clip info ${movie_stackname}_invflt.mrc
  #
  echo "# IMAGE-IMPORTANT: ${movie_stackname}_invflt.mrc <Inverse Fourier-filtered Image (2D, with DW)>" >> LOGS/${scriptname}.results
  #
  if ( ${tempkeep} == "n" ) then
    \rm -f SCRATCH/TMPvalue1_fft.mrc
    \rm -f SCRATCH/TMPvalue2_fft.mrc
    \rm -f SCRATCH/TMPvalue2_2_fft.mrc
    \rm -f SCRATCH/TMPvalue1_msk_fft.mrc
    \rm -f SCRATCH/TMPvalue1_msk2_fft.mrc
  endif
  #
echo "<<@progress: 100>>"
#
##########################################################################
${proc_2dx}/linblock "${scriptname} - normal end."
##########################################################################
#
